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A systematic perturbation theory is developed to describe the magnetic field-induced subdomi- 
nant s- and da;y-wave order parameters in the mixed state of a d^2_y2-wicve superconductor, enabling 
us to obtain, within weak-coupling BCS theory, analytic results for the free energy of a d-wave su- 
perconductor in an applied magnetic field Hci^H <^ Hc2 from down to very low temperatures. 
Known results for a single isolated vortex in the Ginzburg-Landau regime are recovered, and the 
behavior at low temperatures for the subdominant component is shown to be qualitatively difi'erent. 
In the case of subdominant d^^y pair component, superfiuid velocity gradients and an orbital Zeeman 
effect are shown to compete in determining the vortex state, but for realistic field strengths the latter 
appears to be irrelevant. On this basis, we argue that recent predictions of a low-temperature phase 
transition in connection with recent thermal conductivity measurements are unlikely to be correct. 

FAGS Numbers: 74.25.Nf., 74.20.Fg 



I. INTRODUCTION 



Vortices in classic superconductors involve a winding of 27r of the order parameter phase 0(r) around the vortex Hne. 
Since the order parameter = lAkle**^ in a d^2_y2-wave superconductor Hke the high-T^ cuprates is also simply 
a complex scalar with single global phase </>, it was initially expected that the vortex state in the cuprates might 
be structurally quite similar to the textbook case. One remarkable difference was pointed out by Joynt in the 
d-wave case, if subdominant pair potential components of different symmetry exist, corresponding order parameter 
components can be induced at Tc by any probe which couples to gradients of the order parameter, even if the 
subdominant zero field "bare" critical temperatures are very small or zero. Thus the structure of a d-wave vortex 
generically involves admixtures of different symmetry order parameters. Secondly, as noted by Volovik ||^, the 
traditional roles played by extended and localized quasiparticle states in classic superconductors are reversed in the 
d-wave case. In classic superconductors at low temperatures, Caroli-de Gennes-Matricon bound states in the vortex 
core dominate the electronic density of states because extended states are fully gapped and therefore unoccupied. 
In the d-wave case, the existence of order parameter nodes inhibits the formation of the bound states (their very 
existence is questionable and populates the extended ones, which are found to dominate thermodynamics at 

low temperatures and magnetic fields. 

Several authors have attacked the d-wave mixed state structure problem in recent years, armed with these 

ideas. Early studies focussed on an isolated d-wave vortex, allowing for an induced s-wave order parameter component 
and solving the Bogoliubov-de Gennes equations on a lattice. |^ Within the Ginzburg-Landau (GL) theory, similar 
results were obtained. The s component was shown to have opposite winding number to that of the parent 

dx2_y2 in the core regions. Far away from the core center, it decays as l/r^, and its winding number becomes 3, P-p^, 
implying that there are four extra vortices in the s-field at large distances from the main core. These results 

were confirmed by numerical solutions of the Eilenberger equations by Ichioka et al. |ll| The possibility of an induced 
dxy component was also allowed for in Ref. which concluded that the structure in this case was similar to that 
of the (induced) s-wave case, except that the induced order parameter at large distances was found to decay more 
rapidly, roughly as and have opposite winding number 5. Koyama and Tachiki | |T^ pointed out, however, that if 
the calculation is done in a gauge invariant manner there is an additional term in the free energy not found by Ichioka 
et al., involving a Zeeman coupling of the field to an intrinsic orbital magnetic moment in a state with structure close 
to a uniform dx2_y2 -\- id^y. They furthermore showed, within a GL framework with coefhcients determined by BCS 
weak-coupling theory, that this term is proportional to the particle-hole asymmetry of the normal metal from which 
the superconductor condenses, and dominates sufficiently far from the vortex core. 

Interest in the possibility of order-parameter mixing in the vortex state was heightened by the experimental ob- 
servation Jl^~p^ of a plateau in the thermal conductivity as a function of the magnetic field H when H is above 
some critical value H*. Krishana et al., in particular, speculated that their observation of a sharp kink at H* 
might be explained by the sudden onset of an out-of-phase dxy component at this critical field. The new high-field 
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state was proposed to be fully gapped, with vanishing quasiparticle transport. There are several difficulties with this 
explanation, which we discuss below, but theorists were nonetheless persuaded to revisit the problem. 

Laughlin p7| then pointed out, in analogy to the quantum Hall state, the peculiar nature of the time-reversal 
symmetry breaking d^^^y^ + idxy state, which appears to be quite different from possible ground states in the case 
of the d, s-wave mixture. He proposed that the development of a magnetic moment coupling to the magnetic field 
might account for the phase transition, and put forward a T- and _ff-dependent free energy functional driving the 
phase transition, and found a critical field H* ^ similar to experiment. This special free energy functional does 
not have a microscopic basis and it is not clear yet whether there is such a field-induced secondary phase transition. 
To date, no magnetic-field induced transition has been found in relevant numerical studies of the vortex lattice. [pT|J^] 

Laughlin's argument ignored the physics of the core region, thought to be negligible, but Ramakrishnan pointed out 
that quite similar effects are to be expected due to a combination of superfluid velocity Doppler shifts and Andreev 
reflection near the vortex cores. |l8| In higher fields, he proposed that these local d^^-y^^ +idxy patches might overlap, 
causing a transition to a uniform gapped state. This scenario is similar to one proposed by Movshovich et al. in 
the related case of magnetic impurities in a d-wave superconductor in zero field. Finally, Balatsky |20[| has recently 
investigated the effect of the orbital Zeeman term near the upper critical field, and argued that the d-wave state is 
always unstable to a c?^2_j,2 -|- id^y mixture. An unusual collective clapping mode in association with the relative 



phase of dxy to dr^2_y2 was predicted in this superconductor. |21 



Several important questions have not been addressed in the analyses of these issues thus far and have motivated this 
work. First and foremost, we would like to understand whether a phase transition of the type proposed by Krishana 
et al. is possible. Analyses in the GL regime are not applicable, and numerical calculations [ptPj9|,pT[| are not always 
useful to understand competing physical effects. We have therefore developed a systematic calculational approach 
capable of treating carefully the relevant quasiparticle states in the presence of spatially varying superfiow together 
with the relevant subdominant order parameter components on an equal footing in the low-temperature phase. Our 
theory works for H <C Hc2, in which case the vortex core region can be safely neglected. Secondly, we would like to 
understand the structure of the vortex state and the role of the quasiparticles as a preliminary to the yet-unsolved 
problem of quasiparticle transport in applied magnetic field. In addition to being inapplicable at low temperatures, 
the GL calculations on which most of one's intuition for this problem is based are unable to predict magnitudes of 
physical effects since they are based entirely on symmetry considerations. This is particularly important in the case 
of the orbital Zeeman coupling in the dx2_y2, dxy mixing problem. The magnitude of the induced orbital moment is 
a very difficult quantity to estimate properly, as one might a priori deduce by analogy to the intrinsic orbital angular 
momentum problem in the '^He A-phase. In this case, it was found that naive calculations dramatically overestimated 
this effect, and we show here in fact that the orbital Zeeman coupling in the current problem is quite small. 

In this paper, we adopt a semiclassical approach, expanding the BCS free energy in powers of the local superfiuid 
velocity, local order parameter magnitude fluctuations, and their gradients. We thus neglect states possibly localized 
in the vortex core, and other quasiparticle bandstructure effects in a periodic vortex lattice discussed recently by 
several authors [p2|-p5[. Initially Volovik Q and Kopnin and Volovik ||2^ proposed that a semiclassical analysis of 
this type should be valid only down to a scale {Aq/Ef)\/ H/Hc2- Recent numerical work [^,^ indicated, however, 
that the true crossover scale is much smaller for realistic systems with Aq/ Ep <§; 1. Such flne details of the true 
quantum quasiparticle band structure will also be smeared out by impurity effects. We believe, therefore, that our 
neglect of the vortex core and quasiparticle bandstructure will be justifled for cuprate superconductors at low fields 
{H <^ Hc2) and temperatures, and that the current analysis will thus be adequate. 

We begin by presenting in Sec. II the method, which involves a functional integral representation of the BCS 
free energy F, which we then expand in powers of slow superfluid velocity gradients and small subdominant order 
parameter components. Analytical results for F in the GL regime, the low temperature regime, and, for s-wave case, 
an ultralow temperature regime where nonlinear superfiow effects dominate, are given. This allows us in Sec. Ill to 
calculate the order parameter fluctuations directly. We then apply these results to the comparison of structure of a 
single isolated vortex with s or dxy subdominant pairing at various temperatures in Sec. IV, and go on in Sec. V 
to discuss the prospects for observing a low temperature fleld-induced transition of this structure. In Sec. VI we 
discuss existing experiments and make some comments on the various available scenarios. In Appendix A, a detailed 
derivation of the free energy is presented, while Appendix B is devoted to a general calculation of the spontaneous 
magnetization in a dx2-y2 + idxy-vf&ve superconductor. 
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II. FREE ENERGY 



We start from a two-dimensional (2D) phcnomcnological BCS mean-field Hamiltonian in the mixed state: 



Hmf^Y. I ^'^4(1-) ^ ^ -*V - -A(r) " ^^i\c,{v)+ I <fv<fv' A(r, r')4 (r)cj (r') + h.c 



2m 



j (fvcPv'V{Y~v')\b{v,v')\'^, 



where 



(1) 



(2) 



with <i>i characterizing the irreducible representations dx2-y2, s, and dxy, for which the lowest order basis functions 
over a circular Fermi surface are 




i — dr 



(3) 



A(r, r') = V{t — r')6(r,r') is the pairing order parameter with 6(r, r') — (c|(r')c|(r)), and A(r) the magnetic vector 
potential. Throughout the paper h = kg = 1 units are chosen. The magnetic field in the problem is perpendicular 
to the 2D plane, i.e., along the z direction. We assume that Vd < 0, and Vd, Vs and Vxy take such values that in the 
absence of the magnetic field, the superconducting state is of dx2^y2-wa.ve symmetry with A(r, r') = Ao<I>rf(r — r'). 

In the mixed state, A(r,r') takes the form of A(r,r') = e'-'^^^'>A{^,r- r'), where A(R,/9) = Ad(R)$d(p) + 
I?s(R) + Pj;j,(R)$j.y(/9), with A(i(R) = AQ + VdCR.)- The T>i(r) are the magnetic field-induced pairing order parameter 
deviations from their values in zero field. 

The partition function of Hamiltonian ([^) , after making a canonical transformation |p8|| to eliminate the phase field 
of the pairing order parameter (j}{r), is: 



Z = Zgcxp 



J ^Trln 



Mo 



Vi V 



where 



Zo = exp(-Fo/T), 



Fo 



and, in the momentum representation. 



^2^lAdirW , mr)\' , \Vxy{r)\ 



Vd 



(Mo)k,k' - ^ 



-lUJn - Ck 



V 



xy 



(4) 



(5) 



(6) 



-V5(r) • k' -I- -TOw^(r) 



(i>i)k,k' = J d^re'C^'-'^) "- 
(i^2)k,k' = J d^re'C^'^'')-- 



-Vs(r) ■ k - -mw^(r) 



(7) 
(8) 
(9) 



with Wn = (2n + 1)ttT the fermion Matsubara frequency. At = Ao$dk, and 

eA(r) V0(r) 



v.(r) 



2 m 



(10) 
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the supercurrent velocity. In writing down Eq. we have used the relation ^ • Vs(r) = which corresponds to 
the conservation of the supercurrent. Otherwise, mv^ (r)/2 in Eqs. (@) and (|) has no significant effect and will be 
neglected hereafter. 

The free energy resulting from Eq. (^) is 

oo ^ 

F = -TTilnZ = Fo-TTtIuMo+T V — Trp", (11) 

where 

p^i-f5,fi. (12) 



m 

m— 1 



with g the Green-function matrix which, in the momentum representation, is 

iuin+l^k Alt \ 









V 52k 


ff4k / 





3k = ( i = - 1 'x- ■ (13) 



Here, 



W„k = u,l + El = V4 + A|. (14) 

The calculation of the trace of in Eq. (|ll|) can be done by noting that 

(15) 



where k^+i = ki. In most of the bulk region not close to a vortex core, as well as the v^-induced T>i are spatially 
slowly varying functions. Thus we arc allowed to expand Vs{rj) and T>i{rj) in Pk,ki(r) in Eq. ( |l5| ) as power series in 
their derivatives Q, 

v,(r,) ~ e(--^--)-^^v,(r) = v,(r) + [(r^ - r) • Vr]v,(r) + • • • , 

A(r,) ~ e^"-^— ■)-^^P,(r) = P,(r) + [(r^ - r) • Vr]A(r) + ■ • ■ , (16) 

the first few terms of which make main contribution to Trp™. This property enables us to develop a perturbation 
theory to obtain the free energy. The resulting calculation is straightforward but tedious, and is summarized in 
Appendix A. The final result for the free energy with respect to the pairing order parameters F = J d?rf{r), where 
the free energy density /(r), keeping terms up to quadratic in Pj, is 

fir) = /a.W + /,s(r) + f.yiv) + Sf{r), (17) 
|Arf(r)p /■ /^^\ 

Ur) = A2(r)iVo{L;(T,v,)2?;(r) +L';(r,v,)P^'(r) + [c, + /^^(r, v,)][P^(r)]^ + [c, + <(r, v,)] ^'(r)]^, (19) 
f,y{r) = F,{r) + fO^{v) + Al{T)N^{[c^y + ^^^^(T, v,)] [^^^(r)]^ + [c,, + r7^',(r, v,)] [^^'^(r)]^ } , (20) 

F'^^iv) = A^(T)iVo Q°^(r) ^i?(r)P;'^(r), (21) 

^^4(r) = A^(T)iVo {L;,(r, v,)P;^(r) + i;'y(r, v,)P;'^(r)} , (22) 



where ' and " indicate real and imaginary parts, respectively, Vi = Vij t^f^{T\ i = s^d^y, Cg = (— K;iVo) ^ — 2cd, 
Cxy = {—VxyNo)'^ — Cd, with Cd = (— VdiVo)"^, B{v) is the magnetic induction, and tV„k(r), Q^^{T), Li(T, v^) and 
i]i{T, Vs) are defined in Appendix A. /^^ in Eq. (n8) is the free energy density in association with the dominant d^^-y^ 
component in the absence of and "Dxy fs s.nd f^y are the free energy densities for the s and d^y components, 
respectively. Sf involves terms of high orders in v^, "Di, derivatives of Vs, and mixed s and dxy terms. Note in Eqs. 
(ly) and (EQ) the quadratic terms in Di, i = s, dxy can be easily reformulated to coincide with the familiar GL form: 
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(cj + ^^-t^)|X',(r)|2 + !zl-J2^[p-(r)2 + 'D*{r)'^]. H] In addition, derivatives of V, terms are absorbed into the powers 
of Pj terms by partial integration, as indicated in Appendix A. 

With the free energy in Eqs. (p7|)-(p^) we are in a position to investigate the vortex state. We wih first show 
general results for an arbitrary Vs(r) distribution pl| ] in Sec. Ill, and apply the theory to the single vortex case in 
Sec. IV. 



III. ORDER PARAMETERS 



A. d-wave pairing order parameter Ad(r) = Ao + I'd(r) 

In studying the dominant dj.2_y2 component, Vg and Vxy can be set to zero since the Ad-Vg and Ad-Vxy mixing 
terms appear in higher orders of Vg or its derivatives. The gap equation 9/^^/9Arf(r) = produces 

Ad(r) r d^k Arf(r)cos^2y> ^^^^ 



Vd (27r)2 



where Wnkir) is defined in Eq. (A6) in Appendix A. At H = 0, Eq. ( p3| ) reduces to 

1 T/T^Y- f '^'^ cos^ 2ip f d^k cos^2^ tanh(i;k/2r) 

which yields the well-known Tc formula Cd — (l/2)ln(2e'''ci;£)/7rTc), with ajB the BCS cutoff and 7 ~ 0.5772 the Euler 
constant, as well as the asymptotic behaviors of the gap maximum 

. , . r 3.54[T,(r - r,)]V2, T » Ao(r) . 

- \ 2.i4r, -o.39(t/t,)3t„ r<Ao(r) • ^ > 

The presence of vortices depletes the pairing order parameter Aq by Vdir). Sufficiently far from the vortex cores, 
either T or Ao(r) is larger than Vs{r)kF in the whole T regime, and thus we can treat Vs(j)kF and 'Ddi'r) perturbatively. 
It is straightforward to show that the normalized real part of T>d{r) is 

~ 3/4+e^m«2(r)/A^(T) A|5(T) ' regime 

2 (In 2) m-»^(r) T V, ^ ^ ^ l\ ('F\ 

- l-[9C(3)T3-4(l„2).^Tm^,^(r)l/A^(T) A-(r) Ao(T) ' lis^F < i < Ao^J j ^ (26) 

3-[>J, = ±i I cos(e+;f )|3]K(r)fcf ]3/A3(T) A;^(T) Ao(T) ' 



where cf is the Fermi energy and C(3) is the Riemann function. The imaginary part of I'd(r) driven by the derivatives 
of the s uper current is unimportant compared with X'^(r) in the whole temperature region, for its driving term f'^^^ 
in Eq. ( A13) is very small as discussed in Appendix A. 



B. Field-Induced E's(r) 

From Eqs. (0) and (|l|) we obtain the gap equation for the s component, dfs{r)/d'Ds{r) = 0, which gives 

L' (T Vs) - L"{T Vs) 

^^(') " ~w:tM^v ^"^'^ " "2[c.+<'(T,v.)]' ^^^^ 

where L'^{T,Vs), mT,Vs), v's{T,Vs) and ?7"(T,Vs) are defined in Appendix A. 

For T > Vs -kp, we can plug Lg and rjs obtained in Eqs. ( Al5 )-( A18 ) into Eq. ( p7|) to get the scaling functions of 
V'sir) and p;'(r), 

^l'(r) ^ fi^) (29) 
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where 



Gs2{d, z) = (31) 



with his, • • his defined in Eqs. (Al£) and ( A2Cl| ) 



It is easy to see that in the GL regime, 



gl . ^ 1 ^Fm[vliv) - v%{v)] 



It foUows that the ratio between the real and imaginary parts of 2?s 



(34) 



is of order miity in most of the region outside the core where Vs{r) ~ 1/ro and dxVgxiT^) ^ ^/^o with vq the distance 
to the closest vortex core center. 

In the case of • < T < Ao(T), 

7^,iowT. ^ (In 2) T eFm[vlJr) ~ v%ir)] 

C-s i- 2 4TAo " 

l-2(ln2)j;^ eF9.i;.x(r) 
- ^ 2('ln2l^ + 

2 ^ ^1"^^'' Ao ^ 4TA, 

So the ratio between the real and imaginary parts becomes 



c.-l + 2(ln2)^ + ^^I^ Ag 



Ps'°"^(i-) n T m[z;L(r)-^g^(r)] 

p,MowT(j.) ^Ao 9.«s.(r) ■ ^ ^ 

Comparing Eqs. ( p^ ) and (^7|) we find an extra prefactor T/Aq is acquired at low T, indicative of a suppressed real 
part of the s component with decreasing T. This is an interesting observation in the present work, the consequence 
of which on the structure of a single vortex will be discussed in Sec. IV. 

For r < Vg • \ip, which can be achieved either by lowering T in a certain spatial position or approaching the core 
region at a certain T, the nonlinear effects dominate over the thermal effects and the prefactor T/ Aq in Eq. ( ^ ) is 
expected to be replaced by Vskp /Aq. So the real part is negligibly small, and the imaginary part at T = is 

f.„T=0( ^ 1 - So Vsir)kF/Ao ep dxVsxjr) 

^ Cs-l/2 + SeVsir)kF/2A, A^ ' ^''^^ 



where 5*0 = J2i=±i |cos(6' + It is interesting to find from Eqs. ( p8[ ) and ( p}q ) that the factor T/Aq at kp < T 
is replaced by Vgkp /Aq (up to some prefactor) at T < • \ip, reflecting the nonlinear effect due to the Doppler 
energy shift. Its effect on the free energy and the penetration depth in the Meissner state of a d-wave superconductor 
was extensively discussed in the previous papers by the present authors. |^7[] 

C. Field-Induced V^y{r) 

As shown in Eq. (^), there are two terms, F^^ and F4, competing in driving the induced dxy component. The 
first term F'^^ , the so-called orbital Zeeman term, jl^ ] can be rewritten as 

i^^°^(r) = -M(r) • B(r), (39) 
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with M(r) the effective magnetic moment associated with the internal orbital current of Cooper pairs in a dx^-y2 -\-idxy- 
wave superconductor. A general derivation of this spontaneous magnetization is pres ented in Appendix B. As also 
discussed there, M(r) is proportional to the particle-hole asymmetry a (see Eq. (A41)), making the Zeeman orbital 
effects very small for general density of states. It is therefore worth pointing out that analyses of the dr^2_y2-dxy 
mixing problem based solely on a symmetry-based analysis of the orbital moment pc[ | may lead to unrealistic results. 
Besides, it is important to note the very weak temperature dependence of the coefficient Q'-'^ of F*^^, as shown in 
Fig. |ll 

The second term, i^4, is in the order 0{v^) and contains driving terms for 2?"^^ as well as T^'xy This term can be 
significant particularly if the system has small or zero particle-hole asymmetry. Note that in the s-wave case, it is 
irrelevant except for very short length scales of order the core size, due to the nonvanishing leading {0(v1)) term. 
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FIG. 1. Normalized Q'^^ as a functions of T /Tc- Cd = 3 is chosen. 



The orbital Zeeman term has been invoked by Laughlin and Balatsky as the effect driving a putative 
transition to a time reversal symmetry-breaking dx2-y2 + idxy state induced by field. It is therefore particularly 
interesting to investigate the relative importance of F'^^ and within the current approach. For fields H ^ Hd, 
the overlap of vortices leads to nearly homogeneous B{r) ~ H in space and hence F'-'^ /iVxy] can be taken as a 
constant over the bulk. On the other hand, F4/\'Dxy\, which scales as {^/r)'^, is strongly space dependent. It increases 
rapidly when approaching the vortex core, but decays into the bulk. Thus we expect a critical radius r* beyond 
which F'-'^ dominates over F/^, but for r < r* becomes more important and determines the structure, r* can be 
estimated from F4{r*)/F'^^{r*) ~ 1. It turns out to grow with increasing H and decreasing T. In the GL and low T 
regime, we have 



r*{T,H) = 



eH/7nc 

\ T J \eH/mc ) ^ \ ^ J \ J 



1/4 



1/4 



/ \ 1/2 
Xl/4 



1/4 



1/2 



GL 
lowT 



(40) 



with \f the Fermi wavelength, ^ — vf/t^^o the superconducting coherence length, and Rh = \/ c/eH the average 
intervortex distance. We summarize the above estimates in Table I, and display schematically the competition between 
the two driving terms in Fig. ^. Note the "low T" results for F4 obtain only down to temperatures above the local 
superfluid velocity Doppler shift Vskp. For lower T, the perturbation calculation for F^ that we present in Appendix 
A. 3 breaks down because of the zero modes in E\^. Instead one has to get a full expression for the local Doppler shift 
and the derivatives of the supercurrents which is suitable for perturbation expansion after integrating over momenta. 
This is apparently too complicated to be achieved analytically and thus requires a self-consistent numerical work, which 
is beyond the scope of the present paper. However, based on the nonlinear results for the s component and Q'~^^ , we 
do not expect anything qualitatively new in the extremely low T case compared with the Vs • kj? <^T ^ Aq{T) case. 



TABLE I. Order of magnitude of terms driving d^y in free energy outside of core region 



Term 


GL regime (Ao < T) 


Low T {vskp < T < Ao) 










No AoV^y 


iVo {^) AoV^y 



7 



The resulting da;j,-wave order parameter, for ^ < r < r* where ^4 dominates, is 

V' (r) !M£2ll> V" (r) ^MElIfl (41) 



where L'^y{T,Vs), L'^y{T,Vs), <j^(T, v,), and Vxy{T,Vs) are defined in Eqs. ( p^ ), ( pg] ), ( |A4^ ) and ( |A4^ ), and their 



asmptotic behaviors can be found in Table II. In the region r* < r < Rh where F'~'^ dominates, the d^y component 
reads 



shown in Tabel II lead to 



where Q'^^{T) defined in Eq. ( |A4q ). Eq. Q and the asymptotic behaviors of Q^^, ri'^y{T,Vs), and ■q'^yiT.Ws) 



V" (r) / ^"('^"^ " - + ^^r^^^ Ao(T) « T 

U«H^)-^-6(ln2)^][c.,-i + 2(ln2)^ + ^^]-^^ T « Ao(T) ' 



It is interesting to estimate the value of r* in a real material. It follows from Eq. (gg) that for T ~ T^, r*^^ / Rh 
■sJ£,/\f{H/Hc2Y/^. For high-T^ cuprates, ^/Af 10, and thus r*^^ > Rh for fields H > 0.017?c2- In materials 
with larger ^/Ai?, r**^^ becomes order of for smaller fields. Since from Eq. ( pO[ ) r* increases with decreasing 
temperatures, it seems unlikely that the orbital Zeeman free energy plays an important role in determining the local 
order parameter in the vortex state for fields in the Tesla range. 

We would like to make some further remarks on the magnetic field dependence of'D'J.y for r > r* at low temperatures. 
V'J.y shown in Eqs. (|2|) and (^) are obtained from minimizing the free energy density fxyi'r) in Eq. (|2^). This free 
energy density is up to quadratic in V'J.y, which is sufficient for low field and for generic c^y » 1. In the case of larger 
field and/or special case of Cxy close to or smaller than 1, one has to include the free energy density term cubed in 
V'^y This term can be easily found to be {V'J.y)^/3 with the spatial dependence of the coefficient neglected. Thus the 
free energy density up to cubed in V'J.y for r > r* at low temperatures can be written as 

fly = -jBV'^y + coiV'^yf + ^{V'^y)', (44) 

where cq = Cxy + ri'J.y{T,\-s) and 7 — AqNqQ'^^ e/mc. Minimizing f^y with respect to V'^y we immediately get 



V'^y^^cl+jB/2~co. (45) 

It is obvious that there is crossover of the linear- i? dependence of V'J.y for jB /2 ^ c§ to square root of B dependence 
of V"y for 7-B/2 ^ Cq. This interesting behavior is shown in Fig. (^. 

In a periodic vortex lattice, the supercurrent field Vs (r) in the London approximation reads 

KeG 

and thus there will be special symmetry points where ~ 0. At these points, a careful examination shows that the 
coefficients L'^y and L"^ in (Eq. (|2^)) vanish, and that V^y is driven entirely by F^^ , as shown in Eq. (p5|). This 
coincides with the numerical result of Yasui and Kita. 

TABLE II. Asymptotic behaviors of L'^y, L"y, Q*^^, r]'^y, and r]"y (See Appendix A for definitions of Vs-dependent coefficients 
U) 



Term 


GL regime (Ao < T) 


Low T [v.kp < T < Ao) 


r' 




{^u[ + + ^[/^)/Agr 


r" 

-'-'xy 




(-i[/{' + ^c/^')/A;Sr 




-(aAF)(2cd - 1) 


-(a/eF)[In(4a;D/Ao) - 1/2 - 6(ln2)(r/Ao)] 


Vxy 


-(l/8)A^(r)/7r^r^ + eFmvi{v)/Tx'T' 


eFmw^(r)/4AoT 


Vxy 


-(3/8)A^(T)/7r^r^ + eFmvi{v)/-K''T'' 


-1/2 + 2(ln2)r/Ao + eFmD^(r)/4Aor 
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FIG. 2. Schematic comparison of and F4. 



Up to now, we have not included the effect of the d^2_y2 order-parameter suppression Dd on the dxy component, 
since we have shown in Sec. Ill A that Vd is neghgibly smaU in the bulk region, in which we arc primarily interested. 
However, the mechanism for a transition to a dx2^y2 + idxy state proposed by Ramakrishnan [[ISl involves precisely 
the interplay between Vd and the supercurrent near the core, leading potentially to local d^^-y^ + id^y patches which 
then overlap at some critical field. Motivated by this suggestion, we examine the free energy terms including this 
effect within our approximation. These terms, which we refer to as the order parameter suppression terms F'~'^^ , are 
obtained from Eq. (A9-A12) (and have been already neglected in arriving at Eq. (20)). It can easily be seen that the 
Vxy component couples to derivatives of either I?d(r) or of Vs in these terms, indicative of a pure nonlocal effect as 
found by Ramakrishnan. However, for Vs • < T, our surprising finding is that up to leading order, terms including 
derivatives of 2?ci(r) vanish, leaving 



F' 



OPS 



d^l 



d^k Ek 008^293 



2?d(r)[VrXv,(r)],2?i',(r). 



(47) 



Since the [Vr x Vs(r)]2 factor is dominated by the vector potential A(r) rather than V0 part in Vs(r), it may be 
replaced by eB{r)/{mc). Comparing Eqs. (47) and (pT|), we see that F'-'^^ in the core region is of the same order as 
F^^ in the bulk (in particular, it is also proportional to the particle-hole asymmetry of the normal state), and may 
be viewed as the leading correction to F'-'^ , if Ao(r) in F'~'^ is replaced by Arf(r) = Ao(T) -I- I?d(r). It is therefore 
clear from our previous discussion of F'-'^ that F'-'^^ also gives in fact a very small effect even near the vortex core. 




4 6 8 10 

YB/2Co 

FIG. 3. Magnetic field dependence of 'D'J.y for r > r* at low T. 



IV. SINGLE VORTEX 



The results obtained in Sec. Ill enable us to compare the subdominant order parameters, at various temperatures, 
in the presence of vortices characterized by a certain superfluid velocity field (r) . |^ Our purpose in this section 
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is to test these results in the concrete case of a single isolated vortex. 

We use cylindrical coordinates r = (r, 6) with the origin located at the vortex core center. The phase field of the 
order parameter (j>{r) — 6 leading to V4>{r) — —{\/2mr)0. In the spatial regime of interest, the magnetic field is 
roughly homogeneous and thus A(r) = BrO. For r <C Rh, |V0(r)/2rn| » ev4(r) /(mc), so we can neglect A, and 
simply write the supercurrent components and their derivatives as 

sin0 COS0 sin 20 

Vsx^T, > ' O^Vsxij)^-- (48) 

Zmr Zmr Zmr"^ 

except for studying the orbital Zeeman term, in which case ^(r) is no longer negligible because V x A enters. 



High T Low T 




FIG. 4. Single vortex structure with d^2_y2 and s-wave symmetries only in the GL regime ("High T") and for Vskp <C T <C Ao 
("Low T"). Cd — 3 and = 4 are chosen throughout, and distances are given in units of ^. a) and b) Relative phase of s and 
d^2_y2 order parameters with long arrows corresponding to d^2_y2, short to s, and angle between them to relative phase. c)-f) 
Arguments and magnitudes of normalized subdominant s order parameter "Ds- 



A. d^2_y2-s mixing 

In the GL regime, the s-wave subdominant order parameter has been investigated by many authors, [Q-^ with 
substantial agreement. The existing results in this limit are easily shown to be recovered in the present theory. 
Inserting Eq. (E3) into Eqs. (|33) and (^ leads, in the generic case of 1, to 



Pf^(r) ~0.6c;i( M (cos20 + 2isin20) = 0.3cji( ^ ) {3e^'^ - e-^'^) (49) 
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which is consistent with the results in Ref. (7|-||] in the bulk asymptotic region r ^ ^. As a result, the s component is 
of four-fold symmetry, and the relative winding of s and d components is uniform across the whole vortex as shown 
in Fig. ^. 

In the low T case, as discussed in Sec. Ill B, ^^(r) is smaller than P"(r) by a factor of (In 2)T/Ao, and the vortex 
structure at low T is thus expected to be qualitatively different from that in the GL regime. From Eqs. (^), (|3^), 
and (|3^), one finds that 



piowi (j.) _ 23c; 



T 

(ln2)— cos26' + 2i sin 26 

^0 



(50) 



The magnitude and phase of 2?j,™T are shown in Fig. ^. The relative winding of s and d components at low T takes 
place in a very narrow region of real space near antinode directions set by max(T/Ao, Vgkp / ^o)- In Fig. ^, we show 
the T dependence of and V'^ in a spatial position in the bulk. 



0.004 ^ 



0.003 



0.002 



0.001 



0.2 



0.4 



|D"J/A„(T) 
|D'J/A„(T) 



0.6 
T/T„ 



0.8 



FIG. 5. Real and imaginary parts of Va as functions of T/Tc for Vs • Icf < T. r/^ = 10, ^ = 7r/8, Cd = 3, and Cs = 4 are 
chosen. 



B. d^2_y2-dxy mixing 



As discussed in Sec. Ill C, the competition between F'^^ and F4 divides the outside-vortex-core region into two 
rings: the outer ring region r* < r < Rjj is dominated by F^^ resulting in a rigid d^2_y2 + idj.y superconducting 
state with spatially nearly constant V^y obtained in Eqs. (|4^ ) and (^3|); In the inner ring ^ < r < r*, is more 
important and a spatially varying T>xy is expected. Now w e show t he si ng le vo rtex st ructu re in the inner ring region. 
We first insert Eq. (H) into C/(v,) defined in Eqs. (|A33D , ( |A34|) , ( |A2^ ), ([A2S| ), and ( [A37| ) to find that 



Equations (|l|), (|l|) and | 
vortex case. 



C/{(v3)=8C/^(v,) 



-32 [/^(v. 



2e| sin 40 



t/r(v.)^-^^^'(v.).-^ 



(51) 
(52) 



together with the asymptotic behaviors of L^y in Tabel II imply that in the single 



L'^Y'{T)[^] sin40, L'' ^Y" {T) ] cosW 



(53) 



where Y' ~ -7.4 in the GL regime, and -7.lAo(T)/r for T < Aq, while Y" ~ 9.8 in the GL regime, and 25.9Ao(T)/r 
for T <C Aq. Equations (^ij), (53) and the asymptotic behaviors of ri^y in Tabel II yield 



V'xf^iv) 0, 3.7 
V'^^'^ir) ^ -4.9 



A2(T) ^ eprnvUr) 



8(7rT)2 (7rT)2 
3A2(T) , eFmv^r) 



sin 40, 



8 {nTy 



-1 .^X4 



COS 40, 



(54) 
(55) 
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and hence for 



> 1 



-0.61 ic. 



xy 



(56) 



This result coincides with that of Ref. ||ll[ 
At low T, we obtain 



'^xy^'^{-r) - -12.9 



ei?mWg (r) 



4Aor 

1 



Ao 



T \r 



sin 40, 



r epmVg(r) 



n -1 



leading, for c^y ^ 1, to 



-ll ^ 

T 



(6-^4" + 1.76 e'^'') 



T Vr 



cos 4:9, 



(57) 
(58) 

(59) 



In Fig. ^, we show the winding of d^y component and the magnitude of the normalized V^y at high T (GL) and low 
T respectively. The eight-fold symmetry is more obvious at low T. 



Since many current theories postulate a homogeneous 



idxy state without justification, it is interesting 



to consider the size of the spatially varying part of Dxy relative to its homogeneous component. The smallest 
relevant value of the ratio T)^^ /V'^y , where the superscripts F4 and OZ indicate the relevant driving terms in 
Eqs. (^H and (^3|), respectively, is attained at r = Rh in the physically relevant regime Rh < r* . The value is 
{Ep / Aop{H/ Hc2) in the GL regime, and {Ep / /S.qY{H / Hc2){^o/T) at low temperatures; a simple estimate then 
shows that the spatially fluctuating component is always at least an order of magnitude larger than the homogeneous 
component at experimentally relevant temperatures and fields. 



r 

V;0 ■. 



(a) 





FIG. 6. Single vortex structure with d^2„y2 and d^y symmetries only, (a) Relative winding of dxy (short arrows) and ci^2_y2 
(long arrows) components in both the outer and inner rings, (b) Relative phase of d^y to ci^2_y2 component in GL and low T 
regimes, and magnitude of the normalized V^y in (c) high T and (d) low T, respectively, in the inner ring ^ < r < r*. Cd = 3 
and Cxy = 7 are chosen. Distances are given in units of ^. 
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V. IS THERE MAGNETIC FIELD-INDUCED PSEUDO-PHASE TRANSITION? 



The results presented in Sections III and IV suggest that there always exist field-induced s and d^y components in 
a parent dx2-y2-waye superconducting state where the electronic interactions in the s and dxy channels are nonzero. 
These subdominant components are spatially inhomogeneous and, in a general case of small Vg and Vxy compared with 
Vd, are negligibly small as far as any bulk physical quantity is concerned. However, there is also a special situation in 
which either Vs or Vxy is nearly degenerate with Vd, leading to a small Cg or Cxy In this case, the denominator(s) of 
Vi,i = s, dxy may vanish at some critical temperature T*^, a singularity marking a second pseudo-phase transition into 
a dx2-y2 + T>i state with Dj a homogeneous bulk quantity. Investigation of such a possible pseudo-phase transition 
is particularly interesting in association with the experimentally observed thermal-conductivity plateau as mentioned 
in the Introduction. Since in high- Tc cuprates no such phase transition has been reported to be found in the absence 
of magnetic field, we focus on the question whether there can be magnetic field-driven phase transition. This is 
equivalent to searching for a nonzero T*^ at a finite field which vanishes at zero field. 

We first study the s component. From Eqs. (p2|), (|3^), (^5|), and (|36|), we see that the denominator of can never 
be zero, and the singularity may occur in 2?". In the GL regime, the critical temperatures are 
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TLiys) = -^Cs\l/A - 2eFmz;2/A2(T*J]Ao(T;j. (60) 

It follows that T*g{-Vs) < T*^{Q), implying there is no instability at nonzero field leading to a homogeneous subdominant 
V,. At low T, T* becomes 



4(ln2) 



1 - 2c. 



2T*Ao 



(61) 



Again, T*g{-Vs) < T*g{0) is found, meaning that the magnetic field does not favor such a second-order phase transition 
from a dx^-y-^ to dx2-y2 + is state. Accounting for nonlinear effects does not affect this conclusion. 

A similar analysis can be made in the dxy situation. From the asymptotic behaviors of r]'^y and rj'^y listed in Table 
II, we find that there is also a singularity in "D'^y in the GL regime, but its corresponding critical temperatures are 
smaller than that in V'^y. Except for this point, we reach the same conclusion as the s case: No magnetic field-induced 
phase transition is found. 

At this stage, we may examine the prospects of a possible low-T field-induced first-order phase transition into a 
time-reversal symmetry breaking state (i^2_j,2 -|- idxy as suggested by Laughlin JlTt . The free energy functional of the 
dxy pairing order parameter assumed by Laughlin takes the form 



where v = ■y/2Ao/m. The special T dependences of the second and third terms lead to a weak first-order phase 
transition at T ~ 0.52iy^/2eB/c. We may check the justification of this free energy based on the present theory. The 
first and third terms in the right-hand side of Eq. (^2|) are understood to come from the local free energy F^'''"'^ in Eq. 
([A5|). The second one, as the coupling of the unusual magnetization to the field, is naturally related to the orbital 



Zeeman term F^^ in Eq. ( |39| ) with M(r) shown in Eq. (B8). But the assumed temperature dependence crucial for 
the first-order phase transition, is inconsistent with that of F'^^ , which we have shown to be a weak function of T 
for T <C Ao (see Fig. |^) and a linear function of VJ.y in the limit V^y 0. Therefore, the microscopic calculation in 
the present work does not confirm Laughlin's free energy functional which formed the basis of the first-order phase 
transition found in his work. In the present theory, a small homogenous 2?"^ is found for any magnetic field and 
temperature in the superconducting state. It is swamped by a larger (but still much smaller than Aq) V'J.y spatially 
fluctuating component over almost all of the vortex lattice for physically relevant fields. 



VI. CONCLUSIONS 



In this paper, we have formulated a perturbation theory to investigate the magnetic-field induced subdominant order 
parameters of a clean c?-wave superconductor in the presence of the gauge-invariant spatially varying supercurrent 
field in the mixed state. With the assumption of slowly spatially varying supercurrents and their induced s and dxy 
components in the bulk sample, we are able to derive the free energy as power series in the Doppler energy shift, the 
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derivatives of the supercurrents, and the subdominant components. The free energy is vahd from Tc down to very 
low temperatures, enabling us to compare the resulting s and dxy components at low T with the existing results in 
the GL regime. To leading order, the real and imaginary parts of the s component, driven by local m{v1^ — vly) and 
the derivative of the supercurrent dxVsx, respectively, were shown to have very different temperature dependences. 
In the GL regime, both the real and imaginary parts are in the same order. But at low T , the real part acquires an 
extra small prefactor T/Aq and the resulting s winding happens within a very small region near antinodal directions, 
leaving a rigid dx2-y2 +is state over most of the vortex. It is important to note, however, that this structure does not 
imply a gap in the quasiparticle spectrum, due to the small size of the s component compared with the large Doppler 
shifts near the core. 

To leading order, the dxy component is driven by two terms of different symmetries competing over different parts of 
the vortex lattice. The first is the orbital Zeeman term F'~'^ arising from the coupling of the spontaneous magnetization 
to the magnetic field. Its significance is limited by its small magnitude due to the particle-hole asymmetric effects. 
The second driving term F4 scales as (S./r)'^ for ^ <C r ^ Rh- The crossover scale r* divides the bulk region into 
inner and outer regions dominated by the two distinct physical effects. In the inner region, F4 determines the vortex 
structure, and the dxy component has eight-fold symmetry and its relative phase winds 4 times that of the dx2^y2 
component. In the outer region, F'~'^ is more important, leading to a rigid dx2-y2 + idxy superconducting state, where 
the dxy component is spatially nearly homogeneous and a weak function of temperature. We have shown that the 
crossover scale is of order r* = £,{Rh / ^fY^^ in the GL regime, and becomes (Aq/T)^/* times larger at low T. Our 
best estimate for the high-Tc cuprates suggests that r* is greater than the intervortex separation Rh for fields above 
0.01Hc2, so that it appears that vortex lattice structure in fields of order Tesla is governed by F4 and the orbital 
Zeeman effect is irrelevant. The relative phase between the dx2^y2 and dxy components, as well as the magnitude of 
the induced dxy, are therefore strongly space dependent everywhere in the sample. 

Our results have implications for several scenarios which have been proposed to create a state without quasiparticle 
excitations at low temperatures via the creation of a finite out-of-phase subdominant pair component. No such bulk 
state is found for generic values of the pair potentials Vs and V^. It remains possible that induced core dxy patches 
overlap with increasing field, as proposed by Ramakrishnan [Q, leading to a gapped state at high fields beyond the 
scope of our analysis. However, analyzing the interplay between order parameter suppression I'd(r) around the vortex 
cores and Vs(r), we found it to be a very small effect on the induced dxy component; it therefore seems unlikely that 
the field scale where such a transition may occur can be significantly less than Hc2- 

Finally, we searched for a possible second order magnetic field-induced phase transition for special values of the 
coupling constants, and reached a negative conclusion. No such phase transition into a bulk dx2^y2 +s or a dx2^y2 +dxy 
state is found unless the transition has already taken place in the absence of the field, which is apparently not the case 
in the high-Tc cuprates. Examining the Laughlin free energy driving a first-order phase transition, ( [iTt we found that 
the crucial field-dependent term has an assumed temperature dependence inconsistent with the BCS theory; thus this 
phase transition picture is not supported in the present work, consistent with numerical results of Yasui and Kita. 

The apparent phase transition observed by Krishana et al. |lj] has not been clearly reproduced by other groups, 
and the possibility exists that the effect is due to inhomogeneously trapped flux. It is still interesting, however, to 
ask what kinds of intrinsic phase transitions might be possible in a d-wave superconductor. We have shown that it is 
unlikely that any bulk phase transition can be induced by a magnetic field, at least in the low-field regime where our 
approach is valid. Since our model neglects the vortex core regions, it is conceivable that vortex core transitions such 
as those observed in the "^He system might still be relevant. However, since the number of bound states in the core 
region is small and possibly zero for the present case, it seems unlikely that such a transition would have an important 
effect on the quasiparticles responsible for heat transport at low T. A final possibility, currently under investigation 
within the present framework, is that transitions occur in the vortex lattice structure as a function of field. 
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APPENDIX A: DERIVATION OF FREE ENERGY 



In this Appendix we derive the free energy F in Eqs. (17- p2) w hich is vahd for space region where Vg varies slowly. 
We begin with rewriting Trp™ in Eq. ( pl| ) according to Eqs. (15) and (|l6|), 



Trp" ~ (Trp")(^°'=) + (Trp")^'^"'' , 



(Al) 



where (Trp™)('°'^^ corresponds to Vs(rj) and T){rj) for all j — 1, • ■ • ,m — 1 taking the local values v,j(r) and T){r), 



respectively, and becomes 

and (Trp™)*^'^°'^^ contains derivatives of "Vsi^j) and/or V{rj). Inserting Eq. (|Al[ ) into Eq. ([TTI ) leads to 

where F^^°'^'> is just the local free energy obtained in the semiclassical approximation, 
^^(1-' -Fo-TTrlnMo + T 1 (Trp" - / ^'rT^ / ^In hW?„k(r) - ryk(r)], 



(2^ 



i?2(r)+7?k(r)- 



)]/T 



with 



W^„k(r) = -[ic^„ + v,(r) ■ kpf + El{r), 
'^kW^ 5] {[A(r)P*(r) + /i.c.]$dk<i>.k + |A(r)P$2k} 



(A2) 
(A3) 

(A4) 
(A5) 

(A6) 
(A7) 



Expanding F^^°'^^ in power series in Vi will give linear-D^ term as well as quadratic of Vi terms. The prefactors are 
not universal in the whole temperature regimes, as will be discussed and shown below. F'^'^'^''^ in Eq. ( [A^ ) reads 



^(dcr) _ rp _ rjj- p™'\i<i°'^) — y^^(j)^ 



(A8) 



m— 1 



where we expanded (Trp™)''^'''') as power series in the jth derivatives of or Vi with respect to r: (Trp'")^^'^'') = 
^j{TTp"^YJ\ ^(dor) j-ggects nonlocal couplings of the subdominant order parameters to the supercurrent fields. We 

examine the formal leading term, F'-^^ which includs VrWs(r) or Vr2?i(r), in Eq. (Rq), 



oc oo 



m—1 
m-2 



'It p JZl /* 



Pm2 = « 



Pm3 = « 



/ ^'rTr{[VrPk,k(r)] •C^2-''(r).[VkPq,k(r)] -K.kW} 
(m- 1)1^1 d^rTr{[(Vr-Vq)pk,q(r)]q=k -Ck'W}- 



(A9) 

(AlO) 

(All) 
(A12) 



can be expanded as power series in Vi. Note in terms including Vr2?i(r) the derivative can be transfered to that 



of Vs by partial integral. The resulting linear- 1?^ term in F*^^^ is F^^ = , „ F^^\ where 



15 



m — 1 



'^'^^E / 72^*»''{^"-[v.(r) . kj.]} . J2 [-v.(r) • kj.]™-^ J] ^^Tr {(Vk^k) • sl^ ' A(r) • 5^"'^"'} 

n ^ ' m=2 i/=l 

= y d2rPnr)T^ J ^ci>,k2Ak[Vkek-Vr][v,(r).k;^]Z„k(v.) 
-y d2rPf(r)r^ y ^<i>,k2ek[VkAk-Vr][v,(r).k;^]Z, 



nk(,Vs;, 



(A13) 



with r'i(r) 







V*{v) 



and 



{[V,(r) • fcH2 _ ^r^^Y _ 4^,2 [^^(j.) . fc^j 



{[(v,(r) • kpY - M^„k]2 + VJv,(r) • fcj.]2}^ 



(A14) 



It is easy to see that p'f^ is neghgibly smaU, so Fj^ in Eq. ( p^ ) is simply obtained from f ('°^) in Eq. (A4) by setting 
Vs,T>xy = 0. We argue that and Fs^^ in Eqs. (A4) and ( Al3| ) are sufficient for studying the s component up 

to leading orders. Expanding as power series in Ds and Vskp, we will see that the driving term for D'^ is in 

leading order of m{vg,^ — ''^sy)7 which, in the spatial regime of interest ^ < r < Rh, scales as Fi^^ gives a driving 

term for P" scaling as dxVgx ~ l/r^. Clearly, terms contained in X];>i higher orders of Vs and derivatives 

of Ws- We first show results for s component. The dxy situation is more complicated and will be discussed later. 



1. Scaling expressions for L3(r, Vs) and 77s(r, Va) at T > • 



T'/rr \ AT-lrr^ I "^^^ 2 COS^ 2 (4^;^ - W^k ) .2 2 1 ^ SF^bL -^syJ ^ 



r2 



In this temperature regime, we can expand W„k(r) in Eqs. (A4) and ( A13 ) as power series in Vgkp, X'^(r), and 
T>'l{r) to find exphcit expressions for Ls{T,Vs) and ?7s(T, Vg) in Eq. (p^, 



(A15) 

(A16) 

(A17) 
(A18) 



r"frr ^ I c^^k 4cos^2(^ 



7?;(r,v,)=77;'(T,r) + 2d2/^2^(d), 
where Eq. (||) was used, d = Aa{T)/T, and 



4^^/.24rf), 



2 1 d^fiv) 

^2rj drf ' 



hn{d) — — / / (i93Cos^2(^$ 

Jo Jo 

/i3,(d) = - / dx I dLp{l - 2 cos^ 2lp) $ 

Jo Jo 



h2^{d)^-J dx j difcos'^2^'^^^-^, 



^2 tanh(77/2) 

ik" 



2?7 



/i4i(d) = - / dx I dip<i>fy, 
I" Jo Jo 



1 d^fiv) 
2rj drf 



with 77 = y^2-:rd2"cos22^, /(ry) l/{e^ + 1), and = 1 - 2f{r]) + 2-qdf {-q) / d-q . 



2. L,(r,v,) and 77,(T,v,) at r = 



(A19) 
(A20) 



At extremely low T, the perturbation expansions of the free energy as power series in Vs{'r)kp and 'D'i{r) by 
expanding Wnki^) before integrating over k breaks down due to the existence of the zero modes in E'k- The correct 
approach is to expand the free energy after doing the integrals over k. For simplicity, we only show the T = results. 
It is easy to find that the driving ter m f or is negligibly small compared with that for 2?". We may simply set 
= 0. The local free energy in Eq. (A4) becomes 
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Fr 



Hoc) 



T=0 



cPv / A^^El{v)+^^{v) + 2lSv 



rk(r)e(-yk(r)) 



where lk(r) = \j ^'^{v) + [P^'(r)]2 — Vs(r) • kf- and Q{x) is the Heaviside step function. We obtain F^°^^ 
J (fr[cs + ri'JiT = 0, Vs)][P^'(r)]2, where AF is the local free energy at = 0, and 



2 ' ~" 2Ao 

with = Y.i=±i 1008(6* + ;f)|. 

From Fs^"^ in Eq. (A13) we get the driving term for I?", from which we find 



(A21) 



AF ■ 



(A22) 



il'(T = 0,v,) 



y2 - w2 1 

sx syi 



(A23) 



3. Free energy with respect to d^y component 



We first analyze the driving terms at T > • k^?. The leading-order linear-X'a;^ term from i<"('°^) in Eq. ( A4) turns 
out to be 



Moc) 



dW^{T)N^Q'^{T)%{vs)V', 



where 



Q^(r)^A^o-^TX:/(^sinM^{- 



IQEl 



(A24) 

(A25) 
(A26) 



For ^ < r < Rh, F^""^^ is of order 1/r'^. As for Fx]) in Eq. ( |A13| ), the two terms in the right-hand side are qualitatively 
different. We leave the discussion on the second term for a bit later. The first term is nonzero in the leading order of 
{mv1)dxVsx which scales as l/r'* too. For T > Vskp, this term is 



d'rAl{T)NoQ'l{T)U^{^s)V'' 



where 



■ sin 4ip < - 



6El 



rV 

U['{-Vs) = ep {m [ws^(r) - v'^^y{r)] [dxVsy{r) + dyVsx{r)] + 4mVsx{r)vsy{r)dxVsx{r)} 



(A27) 

(A28) 
(A29) 



A simple analysis shows that there are terms coming from F^^^ and F*^^^ in Eq. ( [A^ ) which are of the same order. 
After some algebra, we find that the terms of order l/r'* in F^^^ are 



P(2) 



d^rAl{T)No {Q[{T)U[{vs) + Q',{Tmvs)}'D',,, 



where 



24el 32elEl 



t^i(vs) = e| {vsxi'r) [dlvsy{v) + id^dyV^xi-r)] + Vsyir) [-dlvsxi'r) + idlv sxi'c)] ] 

U2{-Vs) = e\[dxVsx{Y)][dxVsy{v) + dyV^xi^)], 



(A30) 



(A31) 

(A32) 

(A33) 
(A34) 
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and in F^^^ are 



where 



Fif ~ / d^rA5(T)]V„0;(T)J/?(v,)©;' 



24 



116 

— (r) d^Vsxir) + —dldyVsx{v) 

m TO " TO 



(A35) 

(A36) 
(A37) 



The sum of Fj''"^^ and i^iy^ leads to F4(r) term in Eq. (||), with 

L^,(r,v,) = Q;(r)c/{(v,) + Q^(T)[/^(v,) + g^,(T)[/^(v,), (ass) 

L;'^(r, V,) = Q'/(T)C/r(v,) + Q'i{T)U'^{vs). (A39) 

The second term on the right-hand side of Eq. ( Al3| ) is nonvanishing by noting that [VkAk • Vr][vs(r) ■ kp] = 
Aq sin2(y9[Vr X Vs(r)]2 which picks up the magnetic vector potential A(r) in Vs(r). This term is nothing but F'~^^ in 
Eq. @, with 



g-(T).-^-r5:/|^2..sin^2^Z„.(v.).-^-T5: ^ '^"^ '^^^''^''^ 



(2^)2 



(A40) 



pOZ yg].y different from F4. Discussions about comparative importance of these two terms can be found in 
Sec. III.C. The integrand of Q'~'^ includes odd power in ek, implying nonzero contribution only in a particle-hole 
asymmetric system. In order to proceed further, we adopt the following density of states near the Fermi surface in 
the normal state, 



N{e) ~ A^o + e 



dN{e) 



de 



~ A^o(l + a— ), 

e=0 



(A41) 



to take into account both the particle-hole symmetric and asymmetric contributions. Here a is of order unity, and e 
is typically order of T or Aq, implying that the second term on the right-hand side of Eq. ( A41 ) may be negligible 
when the contribution from the particle-hole symmetric mode does not vanish. 
The quadratic free energy terms in Eq. (^) are extracted from 

F(ioc)^ where 



eFTOt;g^(r) 

y2 



hixy{d), 



(A42) 
(A43) 



with hidxy defined in Eqs. ( [A19| ) and ( [A20| ) 



APPENDIX B: SPONTANEOUS MAGNETIZATION IN THE I>x2-y2 +IDxy STATE 
To proceed with a general derivation of the spontaneous magnetization, we rewrite the free energy in Eq. (|ll|) as 



F = Fo - T ^ Tr In (-G^^) - F ^ Tr In 



1 -G 



Vi 
% 



(Bl) 



where G — — (il/o + D) ^ with D the field-induced off-diagonal pairing order parameter matrix. The orbital Zeeman 
term coming from the linear-^ term is 

^-TY^j J d^rd^rie^C'-'')^— ) [G,,, + G^^k.k.W] [(n - r) • V^] [v,(r) • k^] 

-^E / / rf'i-limk,^k{iVk [Gnk,ki(r) + G22k,ki(r)]-Vr}[v,(r).k;.] 

^'^2^^E / (^limk,^k{(^kx Vk)[Gnk.k,(r) + G22k.k,(r)]}-B(r). (B2) 
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Comparing Eq. (B2) with Eq. ( |39| ) we see that the magnetization is just the expectation value of the magnetic 
moment operator {e/2mc){ih x Vk), 



M(r) 



2mcJ (27r)^ 



Tr(4(ik X Vk)4) 



e 

2mc 



(2^ 



hmki^k(ik X Vk)TrGk,ki(r), 



(B3) 



with djjj. = (ck|,cLj^|). It is easy to check that M = for a pure dx2-y2- (d^y-) wave superconductor, and also 
for a dx2-y2 + is-wave one in the order of linear s. For a dx^-y-i + idxy-wave superconductor with d^y component 
perturbatively small compared to the (i^2_j,2 component, we use the Dyson equation. 







Gk,k, 9wS{k - ki) + 5k 
where is defined in Eq. ([l^). Inserting Eq. (B4) into Eq. (B3) yields 

^^"■^ = 4^0^^ J (^'^"^'^ [I?.,(r)54k + P:y(r)5ik] (k x Vk)g2k 



2mc 



rf^k , 4Aoek sin^ 2(/j 

\2 xyV^J' 



(2^)2 



(B4) 

(B5) 
(B6) 



which, up to the leading order of D^y and v^, is consistent with the second term on the right-hand side of Eq. ( Al3| ). 

A particle-hole asymmetric system is required to obtain a nonvanishing M(r). To understand the physics, we note 
M(r) in Eq. (|B^) is in fact 



2 2 ^ 



El, 



{< ~ vl) 



1 



iU!n - i?k 



(k X Vk)52k- 



(B7) 



with m| = (1 + e/j/£'k)/2 and ^k = (1 ^ f-k/E-k)/'2- measuring the particle and hole populations, respectively. Eq. (87) 
can be interpreted as a magnetic moment contributed from number currents of both particles and holes flowing in the 
opposite directions, which cancel in a particle-hole symmetric system. 

The spontaneous magnetization Ajf(r) as a full expression of V'^y can also be obtained from Eq. (B3). After some 
straightforward algebra we get 



M(r) 



d^k 4Ao(T)2?;'^(r)ek<&^,k 



eT X - 



(B8) 



M(r) in Eq. (BS) is important for the direct comparison with Laughlin's free energy as discussed in Sec. V. 
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